Density functional study of the electronic and vibrational properties of TiOCl 
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We present the phonon spectrum of TiOCl computed using hybrid density functional theory (DFT). A 
complete analysis of the spectrum is performed for the space group Pmmn (high symmetry phase) and the 
space group P2i/m (low symmetry phase), which is the symmetry of the spin-Peierls phase. We show that 
the nonlocal correlations present in the hybrid DFT approach are important for understanding the electron- 
lattice interactions in TiOCl. The computed frequencies compare well with those observed in Raman and 
infrared spectroscopy experiments and we identify the origin of an anomalous phonon observed in Raman 
spectroscopy. The relationship between relevant zone boundary phonons in the high symmetry phase and the 
zone center counterparts in the P2i/m symmetry allow us to speculate about the origin of the spin-Peierls 
phonon. 



PACS numbers: 71.15.Mb,63.20.-e,78.30.-j 

I. INTRODUCTION 

The unconventional properties of the low-dimensional 
spin-1/2 compound TiOCl have been a subject of in- 
tensive debate in recent yeo^i^^^^^^^^^^^^^. 
This system shows two consecutive phase transitions 
at Tc2= 91K and at Tc-^= 60K from a magnetic low- 
dimensional behavior to a spin-Peierls dimerized state^i^ 
through an intermediate, structurally incommensurate, 
phase. While the conventional spin-Peierls phase in a 
spin chain is well understood in terms of the magnetoelas- 
tic coupling, the behavior in TiOCl is more complicated. 
The a6-planes of TiOCl consist of TiO bilayers separated 
by layers of Cl~ ions stacked along the c direction^ . This 
bilayered structure, in conjunction with the localized na- 
ture of the Ti 3(i electrons, has motivated various propos- 
als about the role of interchain frustration^ (within one 
layer and between layers), lattice degrees of freedom and 
correlation e^ect^^^^ at the phase transitions. Initial 
suggestions about the influence of orbital fluctuations^ 
have been settled both experimentally and theoretically 
and the groundstate of the system at temperatures above 
is now known to be characterised by singly occupied 
Ti U:,y orbitals'^^^^^^io. 

The ground state electronic structure of this sys- 
tem in the observed room temperature structures has 
been calculated within ah initio density functional the- 
ory (DFT)ii^i^i^. The role of electronic correlations be- 
yond those included in the local density approximation 
(LDA) have also been studied using the LDA+U^^^'^, 
LDA+DMFT^°^^\ and very recently the LDA+cluster- 
DMFT^^ approximations. In the latter approach, where 
two-site correlations were considered, a good agreement 
between the calculated spectral function and the mea- 
sured photoemission and O K-edge X-ray absorption 



spectra was achieved; for all other methods mentioned 
the computed band gap is signiflcantly below that ob- 
served. These studies provide strong evidence for the 
important role played by electron correlations in this ma- 
terial. 

The interplay between the structural and electronic 
degrees of freedom within the various phases, which is 
particularly relevant in the study of the anomalous spin- 
Peierls transition, has only been partially investigated 
from a theoretical point of view. In a previous study 
by some of the authors^, the Raman active Ag phonon 
modes in TiOCl for the high-symmetry (HS) phase (T > 
) were obtained within an LDA+U formalism and the 
frozen-phonon approach. It was concluded that the in- 
clusion of electron correlation at the LDA+U level sig- 
niflcantly improved the agreement between the computed 
and observed frequencies. In the present work we extend 
these studies to include the P2i/m phase, which has been 
identifled by Palatinus et al}^ as the low-symmetry (LS) 
phase (T < TcJ. We present the electronic structure and 
a complete classiflcation of the phonon spectrum for both 
the HS and the LS structures, and make a detailed com- 
parison with the measured Raman and infra-red spectra. 
We also extend the treatment of electron correlation to 
the hybrid exchange DFT approximation - which includes 
both the local (U) and non-local components of exchange 
and correlation beyond those present in the LDA. The 
much better agreement with the observed spectra com- 
pared to that achieved in previous calculations (within 
the LDA and LDA+U)^ strongly suggests that nonlocal 
correlations are important to describe the effects of lat- 
tice vibrations in TiOCl. Within this theoretical frame- 
work we can also identify for the flrst time the anoma- 
lous phonon observed in Raman spectroscop}*^ which ap- 
pears to violate the Raman selection rules. Moreover, 



2 



by analysing the relationship between the zone bound- 
ary phonons in the b direction (Ti-chain axis) in the HS 
phase and relevant zone center phonons in the P2i/m 
phase, we suggest possible candidates for the spin-Peierls 
phonon mode. 

The paper is organized as follows: in section II we 
present the computational details of our calculations, in 
section III we discuss the electronic structure of TiOCl 
in the HS phase within the various approximations used 
in previous studies and present our results using hybrid 
exchange DFT. Sections IV and V report on the opti- 
mization of the lattice parameters and the calculation of 
phonons for the HS phase. In section VI we compute and 
discuss the phonons in the P2i/m symmetry (LS phase) 
and in the last section we summarize our conclusions. 



II. COMPUTATIONAL DETAILS 

The first-principles calculations presented here have 
been performed using the hybrid exchange density func- 
tional B3LYP^^^^^ as implemented in the CRYSTAL 
packagei^. In CRYSTAL, the crystalline wave functions 
are expanded as a linear combination of atom centered 
Gaussian orbitals LCAO with s, p, (i, or / symmetry. 
The Ti, O and CI atoms are described by a triple valence 
all-electron basis sets: an 86-41 IG** contraction (one 5, 
four 5j9, and two d shells), an 8-41 IG* contraction (one s, 
three sp, and one d shells) and an 86-31 IG* (one 5, four 
sp, and one d shells), respectively. A reciprocal space 
sampling on a Monkhorst-Pack grid^^ of shrinking factor 
equal to 6 is adopted after finding it to be sufficient to 
converge the total energy to within 10~^ eV per unit cell. 
The Gaussian overlap criteria which controls the trunca- 
tion of the Coulomb and exchange series in direct space 
are set to 10"^, 10"^, lO"'^, lO"'^, and 10"^"^ the details 
of which are presented elsewhere Typically linear mix- 
ing of 70% and an Anderson second-order mixing scheme 
was used to guide the convergence of the self-consistent 
field procedure. 

The geometry optimizations are performed using the 
algorithm proposed by Schlegel et al?'^ . The full dynam- 
ical matrix was computed using finite differences of the 
analytic gradients for which the atoms were displaced by 
0.02 A. 



III. ELECTRONIC STRUCTURE 

In Ref. 7, the authors performed an investigation of 
the electronic structure of TiOCl at the experimental cell 
parameters of the HS phase (space group 59 Pmmn) us- 
ing the LDA and LDA+U methods. It was concluded 
that the LDA functional does not reproduce the mea- 
sured electronic structure or the internal coordinates 
of the atoms in the cell, whereas within the LDA+U 
method the Mott-Hubbard insulator nature of TiOCl is 
obtained as well as the internal coordinates. Neverthe- 



less, the LDA+U approximation does not describe prop- 
erly the optical band-gap of 2 eV observed in infrared 
spectroscopy experiments. While the gap between the 
occupied spin up dxy band and the unoccupied spin down 
dxy is computed to be 2.8 eV within the LDA+U model 
using U and J values of 3.3 eV and 1 eV respectively, the 
gap between the occupied dxy and unoccupied dxz-yz of 
the same spin is only 0.3 eV. Since the optical gap is a 
measure of the lowest energy difference between occupied 
and unoccupied electronic states, the LDA+U calcula- 
tions significantly underestimate the measured value. 

Saha-Dasgupta et al. improved on the description of 
on-site correlation by using a single-site LDA-I-DMFT, 
where dynamical fiuctuations beyond the LDA-hU are in- 
troduced using a quantum monte carlo impurity solver—. 
This method provides a better description of the overall 
width of the spectral function, but the calculated band 
gap of 0.5 eV still underestimates very significantly the 
observed gap. LDA+DMFT calculations performed with 
an iterated perturbation theory (ITP) impurity solver- 
yield a bigger energy gap at the cost of introducing a 
mixing of the d-orbitals significantly greater than that 
deduced from the spectroscop}*^. 

Based on the fact that TiOCl is a low-dimensional cor- 
related system and that at low temperature it dimer- 
izes into Ti-Ti singlet pairs along the 6-axis, it has been 
suggested^^ that two-site correlations may be impor- 
tant even at high temperatures and thus that the sin- 
gle site approximation will be insufficient. We note that 
in Raman spectroscopy^, specific heat^ and NMR ex- 
periment^Sl strong evidence is found for the existence 
of low energy fiuctuations of the spin system in a tem- 
perature range that extends well above the spin-Peierls 
transition. The role of inter-site correlation has been 
investigated previously in Ref. 12 by using the cluster 
extension of DMFT in which the correlated impurity is 
taken to be Ti-Ti pairs instead of a sing le site-^o. The 
computed linewidth, the lineshape and the band gap of 
the calculated spectral function are found to be in good 
agreement with angle integrated photoemission data and 
with O K-edge absorption spectra^^. 

In the following, we start from a completely different 
approach and demonstrate that the role of nonlocal corre- 
lations can be also captured in an effective manner within 
the hybrid DFT approach. 

The LDA+U scheme is known to lack any description 
of non-local (off-site) electron exchange and correlation 
beyond the local density approximation (LDA) due to 
the fact that the empirical Hubbard U term added to the 
Hamiltonian improves only on the description of on-site 
Coulomb effects. In order to capture the non-local com- 
ponent of the electronic exchange and correlation (which 
is unsatisfactorily approximated within the LDA) one 
would need to extend the on-site LDA+U treatment to 
inter-site interactions, which therefore would take the 
form of an extended version of the LDA+U method, 
namely the LDA+U+V scheme. Within this scheme, 
the inter-site term V is treated within the Hartree-Fock 
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(mean-field) approximation (as tiie U term is in LDA+U) 
and thus the scheme is sufficient to introduce the Ti-Ti 
pair interaction neglected in the LDA. In hybrid DFT, 
the treatment of off-site exchange and correlation effects 
beyond LDA is already built in due to the inclusion in the 
Hamiltonian of a Fock-exchange term which operates be- 
tween all (on-site and off-site) Kohn-Sham orbitals. The 
effectiveness of this term for describing the on-site corre- 
lation is well documented^^ and there is every reason to 
believe that it will be equally effective in describing the 
non-local exchange and correlation effects beyond LDA 
at a mean- field level. 

Fig. [T] shows the spin polarized band structure of Ti- 
OCl calculated within the B3LYP approximation for the 
experimental HS structure^. We note that the band 
structure does not change significantly when calculated 
at the relaxed HS geometry. We observe that, in con- 
trast with the LDA+U approach, the adoption of the 
B3LYP approximation improves significantly on the de- 
scription of the electronic ground state by opening a band 
gap of about 2.2 eV, in good agreement with the mea- 
sured gap. The band states just below the Fermi level 
are majority spin states of dxy character (Fig. [H left 
panel). The coordinate system has been chosen with z 
in the same direction as a, and with the x and y axes 
rotated by 45^ with respect to h and c. Above the Fermi 
level on the majority spin channel (Fig. [TJ left panel) 
we find, in order of increasing energy, band states with 
dxz-yz^ dxz-\-yz^ dx'2-y2^ character and likewise in the 
minority spin panel (Fig. [TJ right panel), with the only 
difference that the dxy spin down states are unoccupied 
and placed between the dxz-yz and the dxz+yz states. It 
is worth pointing out that the bands in Fig. [T] have been 
computed for a ferromagnetic alignment of Ti spins. In 
the antiferromagnetic case (not shown) the band struc- 
ture of up and down states is symmetric, as in any con- 
ventional antiferromagnet, and the band gap is found to 
be the same as in the ferromagnetic case. The states 
at the boundary of the band gap are of dxy symmetry 
(occupied) and of dxz-yz symmetry (unoccupied). 

The above features on the orbital nature of the states 
around the band gap are in full agreement with previous 
results obtained within LDA+lJ^i^. 



IV. LATTICE STRUCTURE 

In Table [J we present the equilibrium lattice parame- 
ters of the HS phase obtained after a full relaxation of 
the atomic structure of the FM state within the B3LYP 
approximation. For comparison, we also show the re- 
sults obtained with the gradient corrected functional of 
Perdew, Burke and Ernzerhof^^ (PBE) and the experi- 
mental values^. 

While the agreement between the calculated and mea- 
sured a parameter is reasonably good, this is not the 




FIG. 1: Calculated band structure for the high symmetry 
phase of TiOCl at the experimental crystal structure along 
the path r-X-S-Y-r-Z. (r[0,0,0]; X[l/2,0,0]; S[l/2,l/2,0]; 
Y[0, 1.2,0]; Z [0,0,1] ) The majority and minority spin bands 
are shown in the left and right panel, respectively. Ep indi- 
cates the Fermi energy. 



case for c and h. The strong disagreement in the case 
of the c parameter is due to the predominantly van der 
Waals type of binding between bilayers in the c direction 
which is not correctly described in LDA, GGA or B3LYP 
approximations to DFT. Keeping c fixed at the experi- 
mental value and fully optimising the structure does not 
change significantly the values of the a and h parameters. 
We notice that the energy landscape is rather flat when 
varying the length of the c lattice constant due to the 
very weak binding. However, we find no significant ef- 
fect of the c parameter on the assignment of the phonon 
modes in the Raman and infrared spectra in the range of 
frequencies considered in the next section. 

Concerning the parameter h (chain direction), the 
B3LYP functional gives an overestimation of h of about 
6 %. In order to understand this issue, we employed first 
the PBE functional for the full structural optimisation of 
the HS phase and then compared with B3LYP results. 
With the choice of the PBE functional, we find that the 
structure is unstable with respect to Ti interchain dimeri- 
sation, where neighbouring chains shift rigidly towards 
each other, leaving the intra-chain step unchanged and 
producing an alternating interchain Ti-Ti distance. This 
is due to the fact that the electronic structure is found 
to be a fully spin polarised metallic state, due to the 
self-interaction error (SIE) present in the GGA. The or- 
bital polarisation found in SIE compensated functionals 
(LDA+U and B3LYP) is lost within PBE and the orbital 
dxy becomes partially unoccupied in favour of the dxz-yz 
orbital. The latter points towards Ti atoms of neighbour- 
ing chains and therefore favours interchain dimerisation. 

This observation is useful in interpreting the behaviour 
of the B3LYP hybrid functional. We note here that the 
functional of Lee, Yang and Parr^^ (LYP) that describes 
the exchange-correlation term in B3LYP, belongs to the 
family of the PBE functional, and thus we expect it to 
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provide similar results as PBE. However, the inclusion of 
Fock exchange removes the interchain instability found 
within the PBE functional and enables the B3LYP ap- 
proximation to describe reasonably well the geometric 
and electronic structure simultaneously. 

By changing the amount of Fock exchange from the 
amount prescribed for the B3LYP functional, we are able 
to conclude that the overestimation of b is strongly re- 
lated to the percentage of the Fock exchange term in 
the Hamiltonian. Possible explanations of the disagree- 
ment between the experimental and theoretical values 
of b are the following: (i) the exchange contribution in 
B3LYP is overestimated due to the neglect of screening 
leading to an overestimation of the local Pauli pressure 
and thus a stretching of b. (ii) more speculatively, if 
dynamical spin singlet fluctuations are important, they 
could bring the Ti atoms into a precursor dimerised state 
with no long-ranged crystallographic order. As a result, 
the x-ray measurement would probe a thermally averaged 
value of the Ti-Ti distance (=b/2) that is shorter than 
the ground state value. It is not clear at present which 
of these explanations accounts for the overestimation of 
the b parameter. 

In Table [TTl we show the relevant equilibrium atomic 
distances and angles. In general a very good agreement of 
B3LYP results with the experimental values is obtained. 

TABLE I: Equilibrium lattice parameters of the high symme- 
try unit cell. 



a (A) b (A) 


c(A) 


exp. 3.78 3.34 


8.03 


PBE 3.98 3.31 


8.57 


B3LYP 3.81 3.49 


8.69 



TABLE II: Equilibrium atomic distances and angles in the 
high symmetry unit cell. 

Ti-0|a (A) Ti-0|b (A) Ti-Cl (A) Ti-O-TiU 



exp. 


1.97 


2.18 


2.41 


147.7° 


PBE 


2.02 


2.13 


2.43 


160.9° 


B3LYP 


1.98 


2.22 


2.45 


148.6° 



V. PHONONS IN THE PACE GROUP Pmmn 

Within space group 59 {Pmmn)^ identified as the high 
symmetry phase of TiOCl, the F point phonons are clas- 
sified according to symmetry as follows: 

F = 3Ag + 3B2g + ^Bsg + 2Biu + 2B2u + 2^3^ 

Each mode involves atomic displacements along a single 
axis. 



TABLE III: Phonon frequencies in the high symmetry phase 
of TiOCl. The first and second column denote respectively the 
irreducible representation and the multiplicity of the modes. 
The third column indicates the direction of the phonon dis- 
placements. The phonon frequencies are given in the fourth 
3-set columns and we include the reported experimentally ob- 
served frequencies in the fifth 3-set columns. 



RAMAN 



Ir.Rep. 


Mult. 


Displ. 


Freq. 


(experim.) 


A, 


3 


c 


205 342 443 


(203) (365) (430) 


B2g 


3 


a 


116 256 613 




Bsg 


3 


b 


141 287 431 





INFRARED 



Ir.Rep. Mult. 


Displ. 


Freq. 


(experim.) 


Blu 3 


c 


— 328 514 




B2u 3 


b 


— 256 274 


— (177) (294) 


Bsu 3 


a 


— 154 463 


— (104) (438) 



The modes with symmetry Ag, B2g, and Bsg are along 
the c , a and b axis, respectively, and are Raman active. 
The modes with symmetry Biu, B2u^ and B^^u are along 
the c, b and a axis, respectively, and composed of one 
acoustic and two optic modes, the latter being infrared 
active. 

In Table UTTl the calculated frequencies are reported and 
classified according to direction of displacement and ir- 
reducible representation. 

The frequencies of the Ag modes at 205, 342 and 443 
cm~^ compare quite well with the experimental values^= 
as do the values computed within the LDA-fU method 
m Ref.0. No cross polarisation of light was considered in 
and therefore no experimental values for B2g and 
B^g are available. 

Comparison with infrared measurements^^ for the B2u 
and B2,u modes, is quite satisfactory for the high fre- 
quency modes but less so for the low frequency ones. As 
a possible explanation for the discrepancy, we estimated 
the anharmonic corrections to the low frequency modes 
by examining the energy surface along the mode and 
found it to be harmonic to an excellent approximation. 
However, it is worth mentioning that the infrared active 
modes have been measured separately by Grueninger et 
al?^ and a value of 201 cm~^ was found for the lower fre- 
quency B2u mode, which agrees significantly better with 
the calculated frequency. 



VI. PHONONS IN THE SPACE GROUP P2i/m 

A. Phonon spectrum 

Within space group 11 (P2i/m), which has been iden- 
tified by Palatinus et al}^ as the spin-Peierls (LS) sym- 
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TABLE IV: Phonon frequencies in the low symmetry phase. 
The first four columns have the same assignment as in Table 
mil In the fifth column we include, where appropriate, the 
HS irreducible representation counterpart. 



RAMAN 



Ir.Rep. 


Mult. 


Displ. 


Freq. 


SG59 


A, 


3 


c 


203 342 443 


A, 




3 


h 


140 289 431 


B3g 




6 


be 


175 212 285 323 382 490 






3 


a 


118 257 616 


B2g 




3 


a 


122 199 540 




INFRARED 


Ir.Rep. 


Mult. 


Displ. 


Freq. 


SG59 


Bu 


3 


c 


— 329 516 


Biu 




3 


h 


— 257 276 


B2u 




6 


he 


175 212 285 323 382 490 






3 


a 


— 159 464 


Bsu 




3 


a 


122 199 540 





metry of TiOCl, the F phonons are classified as follows: 
F = UAg + lOBu + QBg + bAu 

where one Au and two Bu modes are acoustic. 

In this space group, the A^ and modes are Raman 
active while the A^ and B^ modes are infrared active. 
As a general feature, modes in this LS phase are either 
displacements within the be plane or purely along the a 
axis. In particular, the A^ and B^ modes are along the 
a axis and the A^, and B^ are within the be plane. 

We note that starting from the experimental LS struc- 
turei^ within the space group 11 and fully optimising the 
structure results in a structure very close to the HS one 
with a unit cell doubled along the 6-axis. This implies 
that the optimised structure does not show any signature 
of a spin-Peierls distortion. This is due to the fact that 
in the spin-Peierls mechanism the loss of elastic energy is 
overcompensated by an energy gain in spin-singlet fluctu- 
ations. The latter are not properly described within the 
current formulation of spin-DFT and therefore no mag- 
netoelastic structural distortion can be expected for any 
of the exchange-correlation approximations currently in 
use. A consequence of this is that our zone center phonon 
mode calculations for the LS phase provide a good de- 
scription of some relevant zone boundary phonons of the 
high symmetry phase, namely those phonons that are 
associated with a doubling of the HS unit cell along the 
b axis and which are directly related to the spin-Peierls 
instability. 

In Table IIVI the calculated frequencies and the di- 
rections of the corresponding displacement eigenvectors 
(within the reference system a, 6, c of the orthorhombic 
unit cell) are shown. Due to the absence of distortional 
features in the structurally optimised LS structure, the 



18 modes of the HS phase are seen to evolve unchanged 
into the LS phase and the additional 18 modes, originat- 
ing from the doubling of the unit cell, have degenerate 
features which will be discussed below. 

In Table IIVI we choose to decompose the LS modes 
according to their HS irreducible representation counter- 
parts, where meaningful. It is important to note, how- 
ever, that the spin-Peierls distortion does not modify the 
Raman active A^ modes' frequencies of the HS phase 
which are seen to evolve unchanged in the LS Raman 
spectra^. We can therefore assume that the effect of 
lattice distortion is not crucial as far as the frequency 
calculation is concerned. 

A peculiar feature of the two panels of Table ITVl is that 
the strictly LS Raman active modes are degenerate with 
their infrared counterparts i.e. the six A^ with the six B^^ 
modes (third row in the upper and lower panel of Table 
[IV|) and the three B^ with the three A^ modes (fifth row 
in the upper and lower panel of Table IIVP . The six A^ 
and six B^ modes are illustrated in Fig. [2] in order of 
increasing frequency value. 

The degenerate character of these modes is due to the 
perfect triangular geometry of bonds between Ti atoms 
belonging to neighbouring chains that run along the b 
axis. This geometry makes the interchain dimerisation 
patterns shown in Fig. [3] degenerate (this feature should 
be observed in some zone boundary phonons of the HS 
phase). In fact, if we were to consider the case of ladder- 
type arrangement of chains, every Ti atom would have 
only one nearest neighbour in the neighbouring chain, 
and the two dimerisations shown in Fig. [3] would involve 
different lattice strains and therefore would be non de- 
generate. In the present case, where one chain is rigidly 
shifted by half of the lattice parameter b with respect 
to the other, every Ti acquires two equivalent nearest 
neighbours along a and therefore the in-phase and out- 
of-phase displacements shown in Fig. [3] become equiva- 
lent. Note that this degeneracy would be split within 
a model which incorporates spin-singlet fluctuations, as 
pointed out above, although not significantly. However, 
in the non-dimerised phase (T > Td) such a degeneracy 
is restored. 

All the six A^^ and six B^^ modes which originate in 
the LS phase represent possible dimerizations of the Ti 
chains within the be plane. The two lowest energy modes 
(at 175 and at 212 cm~^) keep the Ti-0 distance along a 
(see Table ^ close to its optimised value, while the Ti- 
O distance along b oscillates, in a bond bending fashion 
for the mode at 175 cm~^ and as bond stretching mode 
for the one at 212 cm~^. In fact, the elastic constant 
associated with the Ti-0 bond along a is larger than the 
one along b due to the shorter bond length. The Ti-0 1 a 
bond sets the interchain elastic coupling and, due to its 
presence, the LS modes in the be plane are not merely 
displacements along the b and e axes independently, as 
one would expect if the Ti chains were not coupled. The 
higher frequency modes have larger O components and do 
not generally preserve the Ti-0|a equilibrium distance. 
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175 cm"^ 212 cm"^ 285 cm"^ 




FIG. 2: (Colour online) Eigenvectors of the six and six Bu modes originating in the low symmetry phase, grouped according 
to their degeneracy. Yellow spheres are Ti atoms, red oxygen and blue chlorine. The system of reference is indicated for the 
first mode on the top-left pannel and remains unchanged for the remaining modes. Length and direction of the arrows are 
directly proportional to the strength and vibrational direction of the modes. 
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FIG. 3: (Colour online) Schematic representation of the two 
possible Ti-Ti bond patterns within the P2i/m symmetry 
group for TiOCl. Black and yellow symbols denote the Ti 
atoms on two consecutive layers of a single bilayer. 



as inferable from the length and direction of the arrows 
representing the eigenvectors' components in Fig [2l 

As a general comment on the six and six B^^ LS 
modes, it is interesting to note that they show a low di- 
mensional character. They have a two-dimensional con- 
figurational space {be plane) and due to the strong in- 
terchain coupling some of them show an effective one- 
dimensional oscillation. This can be seen, for example, 
for the modes at 175 and at 212 cm~^ in Fig. [2] where 
the arrows on the Ti and O atoms have same length 
and direction, suggesting that the Ti-O-Ti-0 plaque- 
ttes may be considered as an effective single vibrational 
block. The two modes represent the bending and the 
stretching of the bond between two neighbouring plaque- 
ttes. These features evidence the possible role of low- 
dimensional structural fluctuations as we shall discuss in 
subsection IVI Ci 



B. Phonon assignment and Raman spectra 
anomalies 



We note that the doubly degenerate lowest mode {Ag, 
Bu) at 175 cm~^. i.e., the zone center phonon of the 
LS phase and a zone boundary phonon of the HS phase, 
(Fig. [2] first column, first and second raw panels) remark- 
ably resemble atom displacements along the spin-Peierls 
distortion reported in Ref. 8 and Ref. [27| (compare Fig. 
1(b) of Ref. 8 and Fig. 3 of Ref. f27l with Fig. [2] first 
column, second raw panel). The slight discrepancy in 
the eigenvectors' directions points towards an overesti- 
mation, in our approach of the elastic constant stiffness 
of Ti-Ti bond along 6, which is consistent with the b 
length overestimation in section HVl 

This frequency at 175 cm~^ is close to the frequency 
of an anomalous mode observed in the Raman spectrum 
of Ref. 2 (named a) which consists of a broad peak cen- 



tered around 160 cm~^ at room temperature. As the 
temperature decreases, the peak position shifts towards 
lower frequency values and its broadening reduces until 
it becomes indistiguishible from the Raman A^ mode al- 
lowed in the low symmetry phase. Due to its Lorentzian 
lineshape, this feature is attributed to low dimensional 
energy fluctuations of the spin system which manifest 
themselves in the form of quasielastic scattering of the 
Raman light and decrease as the temperature approaches 
the spin-Peierls transition point^. It is worth pointing out 
that this broad feature is detected when the Raman light 
is polarised along the 6-axis (along which the spin chains 
lie) while it is totally absent when the light is polarised 
along the a-axis^. 

Recently, the evolution with temperature of the Ra- 
man active A^ modes has been studied by Fausti et alM. 
At low temperature, along with the three unchanged HS 
Ag phonons, they find nine additional modes as required 
by group theory analysis (see Table HVl upper pannel). 
Four of these nine modes persist above the spin-Peierls 
transition and become broader losing in intensity as the 
temperature increases. These modes are also visible in 
the Raman spectrum of Ref. 2 at low temperatures and 
the lowest frequency mode coincides with the feature a 
described above. 

Combining the information of the Raman spectrum of 
Ref. 28 (where not only light polarisation bb was used but 
also be) with the information contained in Table IIVI and 
in the phonon displacement eigenvectors (not shown), we 
are able to assign the experimental Raman frequencies to 
the calculated modes (Table |V|). In the following reason- 
ing we assume that the possible splitting of the LS modes' 
degeneracy due to quantum spin effects does not alter the 
picture qualitatively. We have already pointed out the 
fact that the HS Bsg modes remain almost identical to 
their LS (LS Ag) counterpart, as shown in Tables Hill and 
IIVI At high temperature, the Bsg modes are active with 
light polarisation be and therefore we expect their iden- 
tical LS counterparts to have a similar activity. In the 
spectrum of Ref. l28|, the two modes at 178.5 and 524.3 
cm~^ are predominantly active with light polarisation 6c, 
but very weakly active with polarisation bb (this is con- 
firmed in the spectrum of Ref. 2). We therefore assign 
these two modes to the two LS A^ modes with displace- 
ment along b at 140 and 431 cm~^. In order to assign 
the third mode of the same symmetry, we note that the 
modes at 296.5 and 305.3 cm~^ in the spectrum (with 
light polarisation bb) of Ref. 2 are weakly active. Assum- 
ing a systematic underestimation of the frequencies due 
to shortcomings of the theory (which will be discussed 
below) we assign the mode at 305.3 cm~^ to the calcu- 
lated value at 289 cm~^. The remaining six LS A^ modes 
are then assigned straightforwardly in order of increasing 
value (see Table IVl for their assignment) and the three HS 
A^ modes remain unaltered at low temperature (compare 
Table IVl with Table Hill for their assignment). 

By direct comparison, we find that while the three HS 
A^ modes are in good agreement with experiment, the 
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TABLE V: Mode assignment for the twelve Raman active 
modes listed in Table HVl L.S. labels phonons that disappear 
at the SP transition T^ , LP. those that survive up to the 
intermediate phase transition Tc2, H.S. represents the high 
symmetry phase phonons. 



L.S. 



exper. 


145.8 


178.5 


211.5 


305.3 


524.3 


theor. 


212 


140 


285 


289 


431 



LR 

exper. 131.5 296.5 322.6 387.5 
theor. 175 323 382 490 



H.S. 



exper. 


203.5 


365.1 


431 


theor. 


203 


342 


443 



three LS A^ modes along b (Bsg modes in the HS phase) 
are underestimated by 10-20% and the remaining six A^ 
modes of the LS phase are overestimated by 20-30%. 
Possible origin of the above discrepancies could be: 1) 
the lack of a screened exchange term in the Hamiltonian 
which also causes an overestimation of b; 2) the lack of 
spin-singlet dynamical correlations in the B3LYP Hamil- 
tonian which usually renormalise the phononic spectrum 
and, as mentioned in section IV, may be responsible for 
the overestimation of the lattice parameter b. 

Based on the above assignment, we are able to conclude 
that the mode at 175 cm~^ most probably evolves into 
the high temperature a feature of Ref. 2. The remain- 
ing five LS A^ modes show also dimerisation patterns 
which may be tentatively related to the broadening of 
their Raman spectrap^i^^ due to spin fluctuations. Note, 
however, that two of the modes disappear at T^ due to 
a symmetry change in the intermediate phase. It would 
be interesting to perform a phonon calculation within the 
space group 13 (P2/c), suggested in Ref. '28' as the inter- 
mediate phase symmetry, to countercheck the validity of 
our assignment as well as of the suggested space group. 



C. Role of pretransitional fluctuations and the 
spin-Peierls instability 

As seen above, the computed LS phonon calculations 
enable i) the comparison with Raman and infrared spec- 
troscopy measurements performed for the zone center 
modes in the spin-Peierls phase as we showed in the pre- 
vious section and ii) the prediction of the behavior of 
some zone boundary phonons in the HS phase as we will 
discuss in the following. 

From our calculations, the a mode can be interpreted 
as a precursive mode characterised by a strong spin- 
Peierls activity. In fact, this has a very broad lineshape 
and shows a large softening in a temperature range be- 



tween Tci and room temperature^. These features are 
associated with energy fluctuations of the magnetic and 
structural system. Concerning the latter, above the tran- 
sition temperature Tc2 strong incommensurate fluctu- 
ations have been detected in X-ray diffraction experi- 
ments27 and they are seen to persist up to room temper- 
ature. Our phonon calculation supports the lattice fluc- 
tuation scenario with effectively one-dimensional atom 
displacement paths, as shown in subsection IVI Bl for the 
lowest frequency modes. In addition, anomalies in the 
electronic structure are also found in the ARPES spectra 
of the high temperature phase^ which have been given 
little explanation so far. 

Recently, the authors of Ref. |2^ reported the observa- 
tion of a zone boundary acoustic phonon in TiOCl which 
softens at Tc2^^'^^. They identify this mode as the one 
that drives the spin-Peierls transition and estimate its 
bare phonon frequency to be ~ 217 cm~^ . In our 
calculation of the relevant zone boundary modes in the 
HS for the unit cell doubling along b (section VI A) we 
found that the doubly degenerate mode (A^, Bu) at 175 
cm~^ shows features of a spin-Peierls phonon in its vi- 
brational pattern (see Fig. [2]). While we cannot say at 
this point whether this corresponds to an acoustic or to 
an optical zone boundary phonon of the HS phase, its 
energy is close to the bare phonon energy estimated by 
Abel et al^. Moreover, one could then speculate that 
there is a direct relation between the a phonon observed 
in Raman spectroscopy with the interpreted spin-Peierls 
phonon observed in inelastic X-ray experiments. 

A last remark concerns the incommensurate region be- 
tween and Tc^. This region shows incommensurate 
wave vectors shifted from the commensurate wave vec- 
tor (0, ^, 0) by the amount (ibAiJ, ±Ai^, 0). It has been 
argued that the origin of this incommensurability origi- 
nates from competing arrangements of bond patterns in 
TiOCl^ j^^i^^ and more recently in terms of a model of dis- 
commensurations separating commensurate from dimer- 
ized regions^^. The existence of degenerate modes en- 
forces the argument of competing patterns in the incom- 
mensurate region with possible coexisting domains. 



VII. CONCLUSION 

In summary, we have presented within hybrid density 
functional theory a complete description of the phonon 
modes (Raman and infrared active) in TiOCl both in the 
space group Pmmn (high symmetry phase) and in the 
space group P2i/m which is the symmetry of the spin- 
Peierls phase. Comparison with experiment and with 
results obtained from other functionals shows that non- 
local correlations as implemented in the hybrid B3LYP 
approach are important in order to obtain an improved 
description of the modes. This study supports the corre- 
lated nature of TiOCl. 

Our analysis of the zone center high-symmetry modes 
(space group Pmmn) and some relevant zone boundary 
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high- symmetry modes (space group P2i/m) ahows us to 
identify the anomalous a mode observed in Raman scat- 
tering as a probable precursor of the spin-Peierls mode 
strongly coupled to the spin system. Moreover, the tri- 
angular geometry of bonds between Ti atoms belonging 
to neighboring chains running along b leads to the exis- 
tence of six zone boundary modes degenerate with six 
Bu zone boundary modes (for the undistorted system), a 
phenomenon which has not been observed in other known 
spin-Peierls systems and which contributes to the anoma- 
lous behavior of TiOCl. 

Some of the low symmetry modes show effectively one- 
dimensional atomic displacement paths supporting the 
important role of lattice fluctuations in addition to the 



magnetic ones. 
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